Ground-state phase diagram of the quantum Ji — J2 model on the honeycomb lattice 
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We study the ground-state phase diagram of the quantum Ji — J2 model on the honeycomb 
lattice by means of an entangled-plaquette variational ansatz. Values of energy and relevant order 
parameters are computed in the range < J2/J1 < 1- The system displays classical order for 
J2/J1 ^ 0.2 (Neel), and for J2/J1 ^ 0.4 (coUinear). In the intermediate region, the ground-state 
is disordered. Our results show that the reduction of the half-filled Hubbard model to the model 
studied here does not yield accurate predictions. 
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Frustrated quantum antiferromagnets are a subject of 
current intense research, vifhose study is rendered even 
more timely by recent, intriguing proposals of possible 
experimental realizations with ultracold atomsii Frustra- 
tion can arise either from the geometry of the system, or 
from competing interactions, and can lead to the sta- 
bilization of novel, exotic (magnetic and non-magnetic) 
phases of matter. 

A prototypical spin model, featuring interaction- 
induced frustration, is the spin- 1/2 antiferromagnetic 
(AF) Hciscnberg Hamiltonian in presence of next- 
nearest-neighbor (NNN) coupling (usually referred to as 
the J1—J2 model): 




Ji Si ■ Sj + J2 Si ■ Sj 



(1) 



where is a spin- 1/2 operator associated to the ith lat- 
tice site and the first (second) summation runs over NN 
(NNN) sites. Periodic boundary conditions (PBC) are 
assumedi^ 

Model (dl has been extensively studied on the square 
lattice,— and has more recently elicited great interest 
on the honeycomb onC)^"— . For J2 = (i.e., with no 
NNN interaction), the ground state (GS) of ^ features 
AF long-range (Neel) order on both these two bipartite 
latticcSfi^i When J2 > the system is frustrated. Neel or- 
der remains, for sufficiently small values of J2, but as J2 
grows different phases, including disordered ones, become 
energetically competitive. For example, on the square 
lattice the GS displays magnetic long-range order (al- 
beit of different types) for small and large value of J2, 
while the nature of the GS in the intermediate region 
(i.e., J2 ~ 0.5) is still under debate<^ On the honeycomb 
lattice, due to its lower coordination number with respect 
to the square one. the disruptive effect of quantum fluc- 
tuations on magnetic order is enhanced, as confirmed by 
the value of the sublattice magnetization, which in the 
unfrustrated case is approximatively 10% smaller on the 
honeycomb than on the square latticei^. 

Studies of the magnetic phase diagram of the J1—J2 
model on the honeycomb lattice, based on different com- 
putational approaches, have yielded conflicting physi- 
cal scenarios. For example, it is unclear if the GS of 



FIG. 1: (color on line). Schematic representation of the hon- 
eycomb lattice, with the two different definitions of sublattices 
{A and B) corresponding to Neel (left) and coUinear (right) 
orders. 



the system is disordered ("spin liquid") for any value 
of J2, and even the nature of ordered phases remains 
controversial!^""— Exact diagonalization (ED) of ([1]) on 
small lattices (up to 42 sites), yields evidence of AF 
Neel order for J2 < 0.2. For 0.2 < J2 ^ 0.4 the sys- 
tem is predicted to be in a plaquette valence-bond crys- 
talline phase, while, for larger values of J2, a coUinear 
phase (not easy to characterize on small lattices) has been 
suggested^ (see Fig. [1]). The physical picture emerging 
from ED, whose predictive power is clearly affected by 
finite-size limitations, must be assessed by calculations 
performed on systems of much larger size, in order to 
carry out a reliable extrapolation of the relevant physical 
properties to the thermodynamic limit. Series expansions 
yield, for instance, a similar physical scenario, the main 
difference being a GS with spiral magnetic order found 
in the intermediate region.— 

Quantum Monte Garlo (QMG) approaches yield nu- 
merically exact estimates for very large system sizes at 
J2 = 0, but are not applicable in presence of frustration 
due to the well-known "sign" problem. Conversely, vari- 
ational Monte Carlo (VMC) methods arc by definition 
free of such problem, and suffer from no significant size 
limitation. Of course, they are also approximate, their 
accuracy depending on the choice of the variational wave 
function (WF). A recent variational investigation of the 
Ji — J2 Hamiltonian on the honeycomb lattice, predicts 
the disappearance of Neel order at J2 ~ 0.08, a dimerized 
rotational symmetry-breaking phase for J2 ^ 0.3 and a 
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spin-liquid GS in the intermediate region<^ It should be 
noted that, while in studies based on ED the presence 
of order of different kinds is assessed by a computation 
of relevant order parameters, the conclusions of Ref. 
are mainly based on the comparison of energy estimates 
yielded by variational WF's corresponding to differently 
ordered; or disordered GS's. 

Given the qualitative disagreement between ED and 
the only existing variational study, two basic aspects that 
still need be clarified are a) the type of magnetic order (if 
any) of phase(s) that occur as Neel order is suppressed on 
increasing J2, and b) the nature of the GS for J2 ^ 0.4. 
Also of interest is the value of J2 at which the Neel or- 
der vanishes, which according to Ref. is in quantitative 
agreement with that found for the half-filled Hubbard 
mo del 

In this paper, we report results of a variational study 
of the GS phase diagram of ([IJ on the honeycomb lat- 
tice, using the variational family of entangied-plaquette 
states (EPS)ji^ This is a general ansatz, which has been 
shown to yield accurate estimates of GS observables for 
quantum spin lattice Hamiltonians, including frustrated 
oneSfi^ii^ The goal of this work is to determine the phase 
diagram of the model of our interest by direct estima- 
tion of the values of various order parameters, extrapo- 
lated to the thermodynamic limit. Specifically, besides 
GS energies we also compute sublattice magnetizations 
and dimer-dimer correlation functions. We find AF Neel 
order up to J2 ^ 0.2, whereas for J2 > 0.4 the AF order is 
coUinear; for 0.2 J2 ^ 0.4 the GS is instead disordered. 

An exhaustive description of the EPS class of states 
can be found in Refs. UM and We merely wish to 
stress that the EPS WF has been applied to a variety 
of lattice spin Hamiltonians, always providing results of 
accuracy at least comparable to (or better than) that 
afforded by different numerical techniques or variational 
WF's. 

In this study, given the specific magnetic orders (see 
Fig. [IJ that we aim at identifying, the appropriate phase 
factors have been incorporated in the general ansatz. Our 
numerical calculations have been carried out for lattices 
of size as large as = 324 sites. For each lattice size 
and value of J2 , we independently and systematically op- 
timize the EPS WF, by progressively increasing the size 
of our plaquettcs (i.e., the number I of lattice sites that 
a plaquctte comprises). Both the WF optimization and 
the evaluation of various physical quantities are achieved 
via the variational Monte Carlo method. We extrapo- 
late to the thermodynamic limit the estimates obtained 
on finite lattices of different sizes, for a given plaquette 
size, and then, asses the dependence of our results on the 
plaquette size. 

TableUshows the GS energy per site e = E/N for a lat- 
tice comprising N = 100 sites, for different values of J2; 
estimates are shown corresponding to different plaquette 
sizes. Also shown for comparison are the variational es- 
timates reported in Ref. |j by Clark et al. For all values 
of J2 , our energy estimates are lower than those of Ref. 



TABLE I: GS energy per site as a function of J2 for a system of 
size N = 100 with PBC. Data for difFerentplaquette-size I are 
shown. Variational estimates from Ref. |j are also reported 
for comparison. Statistical errors (in parentheses) are on the 



last digit. 


J2 


I = 9 


I = 16 


I = 18 


Ref. 4 


0.0 


-0.54155(7) 


-0.54418(6) 




-0.537901(9) 


0.1 


-0.49049(5) 


-0.49349(6) 




-0.488209(3) 


0.15 


-0.46744(5) 


-0.47087(4) 




-0.469448(2) 


0.2 


-0.44599(4) 


-0.45004(6) 


-0.45145(4) 


-0.450687(4) 


0.3 


-0.40936(6) 


-0.41580(7) 


-0.41751(4) 


-0.41688(1) 


0.4 


-0.40650(5) 


-0.41185(5) 




-0.41125(1) 


0.5 


-0.42432(5) 


-0.42940(6) 




-0.42808(2) 


0.6 


-0.45248(5) 


-0.45712(4) 






0.7 


-0.48519(4) 


-0.48989(5) 




-0.47377(3) 


0.8 


-0.52112(6) 


-0.52610(5) 






0.9 


-0.55928(4) 


-0.56489(5) 






1.0 


-0.59917(5) 


-0.60534(5) 







0. For J2 < 0.1, as well as for large J2, we obtain lower 
variational estimates than those of Ref. already with 
plaquettes of size /=9. Plaquettes of larger size, typically 
/ = 16 (and up to Z = 18 for values of 0.2 < J2 < 0.3), 
arc needed in order to improve on the variational energy 
estimates of Ref. in the remaining cases. This illus- 
trates the main quality of the EPS ansatz, namely that 
it can be systematically improved by increasing the size 
of the plaquette, yielding energy estimates of accuracy, 
typically higher than that afforded by other trial WF's. 

Estimates of the GS energy per site for various val- 
ues of J2 and different plaquette size I are shown in 
Fig. [2] (left part), as a function of the system size. Ex- 
trapolation to the thermodynamic limit of numerical es- 
timates obtained with the same I has been performed 
by assuming the functional scaling form e(A, J2) = 
e(oo, I, J2) + a{l^ J2)N^^/'^. Extrapolated values, there- 
fore, depend on the size I of the plaquette utilized. Since 
by construction expectation values (of any observable) 
must approach the exact GS values in the limit of large 
plaquette size, the question is how to obtain a reliable 
estimate for such a limit, based on those obtained with 
relatively small plaquettes. While in principle calcula- 
tions with ever increasing plaquette size should be car- 
ried out, in order to observe numerical convergence of 
results, currently available computer resources limit the 
largest plaquette size for which this procedure can be im- 
plemented in practice. As it turns out, however, results 
for both the energy and for the relevant magnetic order 
parameters obtained with plaquettes of size up to ^=16, 
lend themselves to simple numerical extrapolations!^ 

An example of this procedure is illustrated in the right 
part of Fig. [51 showing energy estimates extrapolated to 
the thermodynamic limit i.e., e(oo,Z, J2), function 
of the plaquette size. We assume that e(oo, /, J2) can 
be expanded in powers of and fit the data using 

the smallest number of powers. As shown in the fig- 
ure, the values of e(oo,Z, J2) fall on a straight line when 
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FIG. 2: (color on line). Left: GS energy per site of the Ji — J2 
model on the honeycomb lattice with PBC, as a function of 
the system size. Estimates are obtained with A'' plaquettes 
of Z = 4 (boxes), I = 9 (triangles) and I = 16 (circles) sites. 
Each dashed line is a fit to numerical data with the same I 
value. Right: GS energy per site, extrapolated to the ther- 
modynamic limit, as a function of the plaquette-size. Each 
dashed line is a fitting function linear in /^''''■^ (see text). The 
QMC estimate in the J2 = case is also shown for comparison 
(star).i^ 



FIG. 3: (color on line). Left: Neel SSM of the Ji - J2 model 
on the honeycomb lattice with PBC, as a function of the sys- 
tem size. Estimates are obtained with A'' plaquettes of Z = 4 
(boxes), I = 9 (triangles) and I = 16 (circles) sites. Each 
dashed line is a fit to numerical data with the same I value 
(see text). Right: Neel SSM, extrapolated to the thermo- 
dynamic limit, as a function of the plaquette-size. Lines are 
functions bult to infer the I dependence of our data (see text). 
The QMC estimate in the J2 = case is also shown for com- 
parison (star).— 



plotted versus I'^l"^ ^ i.e., an excellent fit to the data is 
obtained with a single powder (dashed lines); in any case, 
the extrapolated value docs not change appreciably on 
including additional terms. For J2 = 0, this procedure 
leads to an estimate in agreement with that obtained by 
QMGi^ii^ (essentially exact in this case). 

The cogent observable, in order to assess the pres- 
ence of magnetic order in the J\—Ji Hamiltonian, is the 
squared sublattice magnetization (SSM), defined as 



m2(iV) = 




where the two summations arc taken over lattice sites 
belonging to different sublattices, as shown in Fig. [T] 
The dependence of the values of the above observable 
on J2, allows one to identify any range of parameters in 
which possibly magnetically disordered phases may exist. 
It is worth restating that in this work we consider only 
two types of magnetic long-range order, schematically 
illustrated in Fig. [U namely Neel and coUinear. 

Figure 13] shows estimates of the Neel SSM extrapolated 
to the thermodynamic limit (left part), based on the func- 
tional form m2(Ar,/, J2) = m2(oo, Z, J2) -h /3(Z, J2)iV"^/^ 



for different plaquette sizes. Just like for the energy, the 
right part of the figure shows the extrapolated values 
m^(oo,/, J2), this time plotted versus The fit to 

the data is performed using a quadratic functional de- 
pendence (dashed lines). 

For J2 = 0, our estimates in the large I limit is 
again in quantitative agreement with QMC one a^'^'^^ . At 
J2 = 0.15 the Neel SSM is finite; we make such a state- 
ment based on the observation that assuming a linear de- 
pendence of m^(cx), Z, J2) on Z~^/^ (the fit is acceptable), 
we obtain a finite extrapolated value [dotted line in Fig. 
[3] (right part)]. On including a term proportional to 
the visible upward bend of the data is better described. 
In turn, the extrapolated value of the order parameter 
increases. 

For J 2 = 0.2 a linear fit gives an unphysical value while 
a quadratic one gives a better description of the data 
leading, within the accuracy of our calculation, to a null 
value of the Neel order parameter. Therefore it is rea- 
sonable to locate the transition point, at which the Neel 
order vanishes, very close to J 2 = 0.2. 

For 0.2 ^ J2 0.4 we find no evidence of magnetic 
order (cither of the Neel or coUincar type), whereas for 
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FIG. 4: (color on line). Same of Fig. [3] for J2 = 0.5. Here 
A and -B sublattices are chosen as in the right part of Fig. 
[T] The dependence of the collinear order parameter on the 
plaquette size is described by using a function which includes 
besides the zero-th order, orders higher than the quadratic (in 
one (dashed line) and only the linear one, excluding the 
point at / = 4 (dotted line). Both descriptions yield, when I 
is saturated, a finite collinear SSM. 
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FIG. 5: (color on line). Finite size scaling of the plaquette 
valence-bond crystal order parameter for two values of J2 at 
which the system is found in a magnetically disordered GS. 
Dashed lines are fits to the numerical data with the same 
J2. This order parameter clearly vanishes, for all the I values 
considered in this work, in the thermodynamic limit. Data 
shown refer to the / = 16 case. 



larger J2, the GS of the system orders coUinearly. 

Figure 2] shows the collinear SSM as a function of iV 
(left part), and its values extrapolated to the thermody- 
namic limit versus I (right part), for J2 = 0.5. The right 
panel shows that the collinear order parameter is finite. 
Collinear order has been found up to J2 = 1.0 (i.e., the 
maximum J2 value considered in this study) 



Our prediction of Neel order persisting up to J2 ~ 
0.2 is in agreement with recent ED works »^ That the 
collinearly ordered phase appears at J2 ^ 0.4 is also con- 
sistent with the ED suggestion)^ As shown by our energy 
estimates (Tab. the collinear phase is favored with re- 
spect to the dimerized rotational-symmetry breaking one 
(described via a resonating valence-bond WF ) proposed 
for J2 > 0.3 in Ref. 0. ED studies also predict a pla- 
quette valence-bond crystal phase in the region in which 
neither Neel nor collinear orders occur. By estimating 
the appropriate order parameter, defined as in Ref. [^, we 
find this kind of order in the region of interest, only in 
systems of small size, i.e., the plaquette structure factor 
is found to vanish in the thermodynamic limit (see Fig. 
[5]). Therefore, our best GS candidate in this parameter 
window, is disordered (spin liquid). A disordered GS has 
been also found by the authors of Ref. Such a phase, 
however appears at a value of the NNN coupling constant 
approximately two times smaller than that found by us. 

Summarizing, the J\ — J2 model on the honeycomb 
lattice has been investigated by means of variational cal- 
culations based on the EPS WF. The phase diagram of 
the model displays three distinct regions: magnetic or- 
der of Neel and collinear type is found for J2 ^ 0.2 and 
J2 ^ 0.4 respectively; in the intermediate region, none 
of the ordered parameters considered here remains finite, 
consistently with a disordered GS scenario. Such a disor- 
dered region has been not revealed by ED studies, in our 
view, due to the small size of the lattices that are accessi- 
ble to them. Interestingly, the value of J2 at which Neel 
order vanishes, according to both ED and our calcula- 
tions, is considerably larger than that (J2 ~ 0.08) found 
in Ref. 0. In Ref. |j, a striking similarity was suggested 
between the phase boundary of the Ji — J2 and that of 
the half-filled Hubbard modeli^ at the Nccl-to-disorder 
transition; the validity of such a suggestion may have to 
be reconsidered, in the light of the results presented here 
which indicate that the simple reduction of the physics 
of the Hubbard model to that of an effective spin-1/2 
system, such as the Ji — J2, might not be achievable. 
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